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Introduction 


This project will address whether the core group of idiopathic autism, often called essential 
autism is neurodevelopmental. We have observed clinically that children in the core population 
have a common facial gestault, which can be described as physically pleasing with wide spaced, 
deep set eyes and a broad forehead. The phenotype is subtle, not generally outside 2SD, so not 
easily catagorized by standard measurements. Study of this facial gestault is considered 
potentially informative since children with complex ASDs are often hyperteloric and facial 
dysmorphogeneis can serve as an index of brain dysmorphogenesis. This project will use 
recently developed 3D shape acquisition technologies and advanced computational techniques to 
define the autism face and determine whether there is a statistically significant facial phenotype. 


Body 
Overview 


We had finished all the proposed tasks proposed in this project (07/01/2008 to 06/30/2010) in the 
statement of work. We had conducted 2D and 3D facial pattern analysis. In the following, we 
will describe each of the finished tasks in details. 


TASK 1: 3DMD® cranial system procurement and installation 


We had obtained $100,000 internal funding from our university and had purchased the 3DMD® 
cranial system (Figure 1) for this project. It is currently installed in the Thompson Center for 
Autism in University of Missouri. 


Figure 1: 3DMD® cranial system (http://www.3dmd.com/). 


TASK 2: _Subjects/Controls recruiting and selection 


We had obtained IRB approval for our project, and had recruited 40 subjects and 72 age- 
matched, typically developing controls. The subjects are recruited from our large clinic 
population in the Thompson Center for Autism in University of Missouri. The age range is 
between 8 to 12 years old. 


TASK 3: __ 3D surface data acquisition and facial feature extraction 


We conducted full 360 degree head/face scan using the 3DMD® cranial system for all the 
subjects and controls recruited for the project. A 3D surface model with both the geometry and 
the co-registered texture image is obtained for each subject (Figure 2). 
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Figure 2: An example of the 3D surface model of the PI with both the geometry (b) and the co- 
registered texture image (a) obtained using the 3DMD® cranial system. 


After the 3D face scan, we use the 3DMD patient® software to obtain 3D coordinate data for a 
suite of anthropometric facial landmarks (Figure 2), following [1]. Figure 3 shows an example of 
using the 3DMD patient® software to calculate the intercanthal distance between landmarks “en- 
en” of Figure 2. 
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Figure 2: Anthropometric landmarks and measurements on face [1]. 
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Figure 3: An example of measuring the distances between anthropometric facial landmarks 
defined in [1] using the 3DMD patient® software. Here the intercanthal distance between 
landmarks “en-en” (Figure 2) is measured. 


TASK 4: Facial phenotype exploration and correlation with clinical observations 


1. 3D facial pattern analysis: 


We had conducted 3D facial pattern analysis for facial phenotype exploration. The research 
result has been published in [4]. We collected 19 anthropometric landmarks (Table 1) from the 
3dMD images from all individuals with the observer blinded to clinical diagnosis. We use a form 
of principal coordinates analysis (PCOORD), which is a clustering method developed 
specifically for landmark data [2]. This analysis enables investigation of how salient features of 
facial anatomy combine on a multivariate axis to define similarities (and dissimilarities) among 
shapes. The goal is to identify the specific combination of morphological variables that 
successfully separate individuals into groups by projecting them onto the multivariate space. The 
computational details of principal coordinates analysis can be found in Krzanowski [3]. 


Phenotype Anthropometric measurement 
Broad forehead ft-ft; ft-ft / tr-g 
Wide nasal bridge en-en / al-al 
Hypertelorism en-en; en-en / t-t; en-en / eu-eu 
Almond-shaped eyes en-ex / ps-pi 
Wide philtrum cph-cph; cph-cph / ch-ch 
Pointed chin go-gn / go-go; sn-gn / go-go 
Outstanding helices sa-sa /t-t 


Table 1 - anthropometric landmarks facial features studied in this project. 


We performed PCOORD analysis on the anthropometric data collected from the faces of the 29 
individuals with ASD for whom we have completed data collection. In order to control for 
effects of size, each individual’s data were scaled by the geometric mean of all possible linear 
distances between landmarks, following. The first two principal axes are illustrated in Figure 4. 
The first and second principal axes account for approximately 19.12% and 17.88% of the 
variation within the study sample, respectively. Note that several individuals on the high 
negative end of the first axis and high positive end of the second axis cluster separately from the 
main group overall (top left of Figure 4a). Suites of linear distances that are highly correlated 
with this location on these axes describe the shape of the mouth and the mouth’s relative location 
within the face as a whole. These individuals display a wider philtrum and wider mouth overall 
relative to the rest of the face as compared to the remainder of the individuals in this sample. 
Also, the mouth is located in a more superior position relative to the upper and midface, such 


that is closer to the nose. As the data have been scaled for size, this is not an effect of age or 
growth. 


(a) (b) 


Figure 4: Facial feature clustering result for ASD based on principal coordinates analysis 
(PCOORD). (a) is the clustering result within the ASD group, (b) is the clustering result by 
combining the ASD group (shown as dot) and the control group (shown as X) together. Only the 
first two principal axes are shown here. See text for more descriptions. 


Though this subgroup of individuals appears to be outside of the range of phenotypes of the 
majority of the individuals in our sample, we cannot know from this analysis whether the 
individuals in the subgroup or the main group are more representative of the typically- 
developing facial phenotype. To answer this question we must include the same data for both 
individuals with ASD and typically-developing individuals. Thus, we performed the same 
PCOORD analysis with the 29 individuals with ASD and included in the analysis data acquired 
from 69 typically-developing individuals. These data were also scaled to control for the effects 
of size. The first two principal axes are illustrated in Figure 4b. First, the majority of individuals 
with ASD cluster together with the typically-developing individuals. However, we see several 
individuals located on the high positive end of the first axis, with this axis accounting for 
approximately 18.65% of the variation in the sample. Again, these data were scaled prior to 
analysis, which indicates that the distribution of individuals is not based on age- or growth- 
related differences; in fact, the individuals included in this cluster are the same individuals that 
clustered separately in the analysis of only the individuals with ASD. Additionally, the suites of 
linear distances that define this subgroup closely match those from the previous analysis, i.e., the 
breadth of the philtrum and mouth, and the more superiorly-located mouth. 


The following are the summary of the study results: 


1) Analyses show that children diagnosed with ASD demonstrate differences in facial 
morphology compared to typically developing children (Figure 5). 


Figure 5: Measures that are significantly different in boys with ASD as compared to typically 
developing children. White lines indicate measures that are larger in children with ASD, while 
purple lines represent measures that are smaller in children with ASD. 


2) There are 2 subgroups of children with ASD that are different from both the other children 
with ASD and typically developing children (Figure 6). These results do not reflect correlations 
with age or head size. 


(b) 


Figure 6: There are two subgroups within ASD. a) Measures distinguishing Subgroup 1. Purple 
lines indicate measures positively correlated with Axis 1, while white lines are negatively 
correlated with Axis 1. b) Measures distinguishing Subgroup 2. Purple lines indicate measures 
positively correlated with Axis 2, while white lines are negatively correlated with Axis 2. 


3) We find statistically significant correlations between subgroup membership and clinical and 
behavioral characteristics (Table 1). 


Table 1: Significant correlations with subgroup membership as compared to the main ASD 
group. 


|) Subgroup 1 | Subgroup 2 


Diagnosis t Autistic Disorder t Asperger's 
Severity High SCQ Moderate SCQ 
IQ Lower Higher 
Attention Higher CBCL Lower CBCL 
LANGUAGE — |] wereennnnnnnnnneene Higher 
Family NiStOry | wwse-en2-222-------- tAnxiety Disorder 
Other t seizures t clumsiness 


2. 2D facial pattern analysis: 


In this study we aim at the following task based on the 2D facial features of each subject: to find 
if there is such MFP group and the homogeneous features shared in this group. Frontal facial 
images of 41 subjects with ASD are used in this study. Sixteen feature points are automatically 
located on each facial image using our newly developed method [5]. To ensure the accuracy of 
the feature point location, all the results are examined by a human operator, and manual 
adjustment can be made if the detected feature point is off its correct location. Procrustes 
analysis is used to align the points across all the subjects. Size difference is removed during 
Procrustes analysis. EDMA is performed which generated a pairwise distance matrix. The upper 
diagonal of this matrix is concatenated into a distance vector. k points has k(k-1)/2 pairwise 
distances. 


Our interest is finding a homogeneous midline facial group (MFP) group in the whole 
population, while the rest of the population is of less interest and may be heterogeneous itself. A 
proper clustering method needs to be chosen in order to suit our application. Traditional 
clustering methods such as k-means may not be a good choice, since we do not have a good 
estimate of the number of clusters due to the lack of hypothesis of the non-MFP group. High 
dimensionality of the feature vectors may also cause problems for many clustering algorithms 
due to the curse of dimensionality. Many dimension reduction methods have been developed 
which include feature transformation and feature selection. Feature transformation represents the 


data by a smaller number of new features, each being a combination of the original features. 
Information from irrelevant dimensions is preserved, which makes these methods ineffective 
when the clusters are masked by irrelevant dimensions. The combinations of the original features 
are also difficult to interpret. Feature selection selects the most relevant dimensions from all 
dimensions, but it may have difficult when clusters exist in different groups of dimensions 
(subspaces). 


Subspace clustering is an extension of feature selection that attempts to find clusters in different 
subspaces. It appears as biclustering in the application of gene expression data, but subspace 
clustering and biclustering are essentially the same. It can cluster the observations and features at 
the same time, thus performing feature selection during clustering process. Subspace clustering 
suits our purpose well, because we are interested in not only looking for an MFP group, but also 
which features can distinguish MFP from non-MFP. 


We use the software BiVisu to perform biclustering. The data are represented in a matrix, where 
each column is a feature dimension, and each row is an observation (subject). 


Classification is performed if we are able to find an MFP group by clustering. The clustered data 
are divided into training and test data. LDA and SVM are used to train the classifier. Since our 
sample size is relatively small and the number of subjects in each group is uncertain, it is 
difficult to divide the data into training and test sets fairly. Instead, we use leave-one-out method 
each time and then calculate the average classification rate. 


Eight subjects have open mouths in their images. Two data matrices are constructed for two 
alternative ways to remove the open mouth effect. The first data matrix contains 33 rows 
(subjects) and 120 columns (distances), with the 8 subjects excluded. The second data matrix 
contains 41 rows and 45 columns, with the 6 feature points on the mouth excluded. Next we will 
show the results of the first data matrix in detail, and later we will summarize the results of the 
other matrix. 


There are several parameters to set. The minimum number of columns of each cluster is set to 
20, the minimum number of the rows is set to 25% of the total number of rows, the maximum 
cluster overlap is set to 20%, and the noise level is set to 0.05. A noise level too high will group 
the entire data matrix into one bicluster, and a noise level too low will find no bicluster. 


Multiplicative model is used and the result gives 8 biclusters. The bicluster with the highest 
correlation is selected as the most significant bicluster. 14 subjects and 34 features are grouped 
into this bicluster. The ratio map of a bicluster is defined as the elementwise ratio of each column 
to a reference column (e.g., 1" column) of this bicluster. Figure 7 shows the ratio maps of the 
selected bicluster. The similar patterns across rows indicate the coherence of this cluster. Note 
that in this experiment there are 33 subjects and 120 features total. 
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Assume that the subjects in the selected bicluster belong to MFP group. Figure 8 shows the 
feature points of each subject, with MFP subjects in red and non-MFP subjects in blue. The 
distances in this bicluster are also shown, where the red lines indicate the average distance is 
longer in MFP group. From this plot we find that the group differences lie mostly in the eyes and 
the area between the nose and the mouth. MFP group have shorter eye width, shorter distance 
between eyes and longer eye height. 
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Figure 7: Ratio map of the bicluster for method 1. 
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Figure 8: Feature points (red: MFP) and the selected pairwise distances (red: MFP>non-MFP) 
for method 1. 


The data are labeled according to the clustering results, and only the 34 features in the selected 
bicluster are used in the classification. The experiment is performed 33 times, each time one 
subject being the test data. The average classification rate of LDA method is 84.8%, and the 
average rate of SVM is 60.6%. 


Now we will show the results on the data set with 41 subjects and 10 feature points (45 features). 
Biclustering results give 4 biclusters, and the selected bicluster has 16 rows and 31 columns 
(Figure 9). Figure 10 is similar to Figure 8 except there are no mouth feature points. Compared 
with Figure 8, the same conclusion about the eye features can be drawn from Figure 10. 
However, the group difference of the eye-nose distance is inconsistent with Figure 8 and the nose 
width difference is not shown in Figure 8. 


The classification is performed in the same way as 4.3. The average classification rate is 88.2% 
for LDA, and 72.3% for SVM. 
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io maps of the biclusters for method 2. 
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Figure 9:.Feature points (red: MFP) and the selected pairwise distances (red: MFP>non-MFP) 
for method 2. 


Summary of Results: 


We obtain two subgroups from the above results. The first is clustered from 33 subjects with 16 
feature points (120 features), and the second is from 41 subjects with 10 feature points (45 
features). By comparing the two groups, we find that they have 11 subjects in common, which 
constitute a majority of each group. Among the rest 5 subjects (16-11=5) in the second group, 3 
of them are not included in the first experiment because of open mouth. 


The two experiments select similar distances among eye feature points, and these distances show 
the same pattern of MFP group (shorter eye distance, longer eye height, etc). A few eye-nose 
distances are also selected in both experiments but the revealed patterns are inconsistent. This 
inconsistency is likely caused by the few non-overlapping subjects of the two MFP groups. 
Those subjects might be the noise in the corresponding clusters, although further experiment is 
needed to confirm it. Nose width is selected in the second experiment but not the first one, and 
several nose-mouth distances are selected in the first experiment only because of the exclusion of 
mouth feature points in the second experiment. 


The subsequent classification is used to test if the selected features have the capability of 
distinguishing MFP from non-MFP. Both experiments show reasonable classification rates. LDA 
performs better than SVM although LDA is simpler in theory. 


Based on what has been done so far, we find that there exists a homogeneous MFP group in the 
ASD population, and this group generally has shorter eye width, shorter eye distance and longer 
eye height compared to the rest of ASD population. Other facial features do not appear to be 
distinguishing in the subgroup, although further experiment is needed to provide more evidence. 
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Key Research Accomplishments 


We had generated a set of precise, highly replicable 3D anthropometric data for core 
ASD children and age-matched, typically developing controls, which has not previously 
been done. 


We had conducted 2D and 3D facial pattern analyses. Results suggest that there are 
subgroups of individuals with ASD that display distinct facial phenotypes. 
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Reportable Outcomes 


We had collected 3D surface model with both the geometry and the co-registered texture 
image is obtained for 40 ASD subjects and 72 age-matched, typically developing 
controls. We will make the data available for the public. 


Results of 3D and 2D facial pattern analyses suggest that there are subgroups of 
individuals with ASD that display distinct facial phenotypes. 
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Conclusion 


We had generated a set of precise, highly replicable 3D anthropometric data for core ASD 
children and age-matched, typically developing controls, which has not previously been done. 
Results of 3D and 2D facial pattern analyses suggest that there are subgroups of individuals with 
ASD that display distinct facial phenotypes. Characterization of the facial phenotype will 
enhance understanding of embryologic forces which can cause autism, and may provide a 
potential prescreening tool to assist early diagnosis. 
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Automatic detailed localization of facial features 


Qing He’, Ye Duan 
Department of Computer Science, University of Missouri, Columbia, MO, USA 65211 
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Abstract 

We propose a complete framework for automatic detailed facial feature localization, for the purpose of computer 
based syndrome classification. Particularly, we are interested in the features of the eyes, the nose and the mouth. 
Face detection is performed followed by the region detection that locates a rough bounding box of each facial 
component, and detailed features are then extracted within each bounding box. Since the feature points lie on the 
shape contours, we start from shape contour extraction, and then detect the feature points from the extracted 
contours. In the eye and nose feature detection, a local minimum map is constructed followed by curve detection and 
fitting. In the mouth feature detection, the lip is first segmented based on the color information and the lip contour is 
extracted by fitting parabola curves. Experimental results show the robustness and accuracy of our methods. 


Keywords: facial feature localization, eyelid, nose boundary, lip contour, Hough transform 


1 introduction 

Automatic localization of facial features is important technique in many applications, such as face recognition, 
human-computer interaction, surveillance, etc. Particularly, computer aided diagnosis is becoming an attracting 
application in recent years. For example, the presence of a specific facial pattern in patients with a genetic syndrome 
indicates phenotypic expression of this disease, and the question has been asked whether a computer can recognize 
and classify dysmorphic faces solely on the basis of 2D images. Dalala and Phadk (2007), Loos et al. (2003), 
Boehringer et al. (2006) performed syndrome classification based on the feature points extracted from facial images, 
but the feature points were manually identified in their works. We aim at developing an automatic method to locate 
facial features for the purpose of computer based diagnosis, but the approach is applicable to other fields with 
similar requirement. 

Various methods for facial feature detection/localization have been reported in literatures. Face detection is usually 
the first step in facial feature detection in order to provide a region for facial features to be detected. Feris et al. 
(2001) proposed a two-level Garbor wavelet network (GWN) to detect eight facial features. In Bhuiyan et al. (2003) 
six facial features are detected. Face detection is performed using skin color segmentation and genetic algorithm. 
Image enhancement is then performed and the facial features are detected through intensity threshold. In Gourier et 
al. (2004), an elliptic facial region is detected based on skin color. A vector of Gaussian derivatives and the linear 
combination of these descriptors act as detection functions for facial features. In Cristinacce et al. (2004), a multi- 
stage approach is used to locate 17 facial features. The face is detected using boosted cascaded classifier (Viola and 
Jones 2004), and the same classifier is trained using facial feature patches. Further refinement is achieved by active 
appearance model (Cootes et al. 2001). Asteriadis et al. (2009) detect eye centers and mouth corners using a 
hierarchical scheme. After face detection and normalization in size, a distance vector field (DVF) is created 
according to the distance from each pixel to the closest edge. The eye and mouth regions are found by matching the 
DVF with a template. Intensity and color information is then used to locate eye centers and mouth corners within the 
detected regions. 

All the above works only look for feature points. However, more details of facial features lie in their shape 
information, such as eyelid and lip contours. While it may not be too tedious to label a few feature points manually, 
it is almost impractical to delineate shapes of facial features by hand. Therefore, an automatic procedure for facial 
feature shape extraction is in demand, which is the main focus of this paper. 

Lip contour extraction has been a branch with extensive studies in the past decade. Since the lip is rich in color, most 
previous works utilized the color information to separate the lip form the background, and then some refined 
methods were used to extract the accurate lip contour. Wang et al. (2004) used three quadratic curves to characterize 
the lip contour, with 16 points on the curves as parameters. The lip region was first separated from the background 
using fuzzy C-means clustering with shape constraints (Wang et al. 2002). Eveno et al. (2004) used several cubic 
curves to model the lip contour. Feature points on the upper lip were first detected using “jumping snake” which is 
robust to the initialization. Hybrid edge (Eveno et al. 2002) was computed from pseudo-hue which intensifies the 
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upper lip edge. Yokogawa et al. (2007) matched a template of four quadratic curves to the lip contour, and the lip 
was first separated from the background by thresholding the hue component. Active contour model was also used in 
several literatures. Wakasugi et al. (2004) proposed an edge detection method based on the separability of two color 
intensity distributions, which was used as the image energy in the active contour model. A level set method without 
re-initialization was proposed in Sohail and Bhattacharya (2007), where only grayscale images were used. 

Much fewer studies exist for eyelid and nose contour extraction. Vezhnevets and Degtiareva (2003) proposed an eye 
contour extraction method. Iris was first detected based on its intensity and circular shape. The points on the upper 
eyelid were then detected by Hough transform. A cubic curve was used to fit the upper eyelid points and a quadratic 
curve was used to fit the lower eyelid. Zheng et al. (2005) used Gabor filter to detect the eye corners and the top and 
bottom points of the iris, which were then used as control points for a spline curve to fit each eyelid. Ding and 
Martinez (2008) performed complete facial feature detection including the shape of the eyes, nose, mouth, eyebrows 
and the lower chin. A classifier was trained for eye corners and mouth corners respectively. The eyelid points were 
then detected and fitted to a cubic curve in a way similar to Vezhnevets and Degtiareva (2003). The mouth and nose 
boundaries were detected using some simple thresholding and the gradient information. Earlier works (Kampmann 
and Zhang 1998; YUILLE et al. 1992) used edge based deformable models to fit the shape contours of facial 
features, but it has been stated in Vezhnevets and Degtiareva (2003) that the edge information is not robust due to 
spurious edges and discontinuous edges. 

We propose a complete framework for automatic detailed facial feature localization (Fig.1). Features of the eyes, the 
nose and the mouth are of interest (Fig.2), and more explanations of these features will be shown later. Face 
detection is performed followed by the region detection that locates a rough bounding box of each facial component, 
and detailed features are then extracted within each bounding box. Generalized Hough transform is applied for both 
eyelid and nose boundary detection, and color information is used for a rough lip region segmentation following by 
a curve fitting. Most of the feature points can be detected from the boundary contours. Our method is different from 
the previous works in that 1) more points on the lower eyelid are detected before the final curve fitting, while 
previous works only used two or three points to fit the lower eyelid contour. 2) it is the first time to use Generalized 
Hough transform for the nose boundary detection, which is more robust than simple thresholding methods in 
previous works. 3) five parabola curves are used to model the lip contour, which have simpler and more flexible 
parameters to adapt to different shape. 

The primary application of our work is computer based diagnosis, so a few assumptions are made on the input image. 
We assume the input is a frontal color face image with neutral expression. The face should take the majority of the 
image area, and should not be occluded by hair, beard, etc. There should be no in-plane rotation or out-of-the-plane 
rotation of the face. Usually these assumptions are easy to satisfy in clinical applications. 

The reminder of this paper is organized as follows. Section 2 describes region detection by parametric template 
matching, section 3~5 describe the detailed methods for the feature extraction of the eyes, the mouth and the nose, 
section 6 shows the results and validation, and section 7 concludes the paper. 


2 face detection and region detection 

2.1 face detection 

Prior to any feature detection, face detection is applied using the method in Nilsson et al. (2007). Briefly, the local 
Successive Mean Quantization Transform features are proposed for illumination and sensor insensitive operation in 
object recognition. A split up Sparse Network of Winnows is presented to speed up the original classifier, and the 
features and classifier are combined for the task of frontal face detection. The input image is first converted to 
grayscale for face detection and the subsequent feature detections, except that the color version is retained for mouth 
feature detection. The image is then cropped to the detected face area as shown in Fig.3(b), which covers all the 
facial components but not the entire face. To facilitate the following template matching, the cropped image is 
rescaled to a standard size, e.g., the size of the template face image. 


2.2 parametric template matching 

Template matching is the process of comparing a standard pattern (template) with an image for the purpose of 
finding the occurrence of the template pattern within the search image. The search is performed by moving the 
template pixel by pixel over the search image while calculating the correlation between the template and the portion 
of the search image at that position to find the position with the highest normalized correlation value. High 
computational cost has been a major difficulty of this method. However, if the search area is small this problem can 
be greatly alleviated. The parametric template (Tanaka et al. 2000) is a linear combination of a set of templates with 
different feature values (scales, rotations, etc). It has many advantages over traditional template matching in that it 
can represent geometrical and non-geometrical changes of an object in the parametric template space. 


Here we briefly introduce the concept of parametric template. A root template ft, is defined as the template directly 
obtained from a template image. A set of vertex templates can be generated from the root template by geometric or 
non-geometric transformation. Let t,,...0,, be M normalized vertex template with the same size, the parametric 
template is defined as: 

t= (wit) wit, | (OSw, <L>"w, =) @ 

When comparing the template and an image patch g (normalized) with the same size of the template, our goal is to 
maximize the normalized correlation O(g,t,,) by varying the parameters W subject to the constraints in (1), 
which can be solved by Lagrange multiplier method. The solution is 
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In order to satisfy the constraints of W in (1), normalization needs to be performed on Was follows. 
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A feature value V, is bound to each vertex template t,. After the parameters Ware obtained, the feature value V of 


the parametric template f,, can be calculated as 


(4) 
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2.3 region detection 

The parametric template method is applied to the region detection. Scale variation of the objects is the major 
concern. Although the face is normalized to a standard size, there is still variation in the size of each individual 
facial component among different people. We take into account the size variation in both x and y directions. 


Therefore, a feature vector V, instead of a scalar is assigned to each vertex template, which consists of the scales in 


x and y directions relative to the root template. 

The root template for each facial component is shown in Fig.3(c). Each template is constructed such that some 
margin is included but no part of other component is included. Five vertex templates are generated from each root 
template, which are scaled in x direction only, y direction only, and both x, y directions. In order to compute the 
linear combination of the templates according to (1), the root template and its derived templates are first aligned by 
their centers, and all the templates are padded to the size of the largest vertex template so that they have the same 
size. The resulting parametric template has the same size as the largest vertex template. 

In the search process, we utilize the prior knowledge about the approximate positions of facial components. Thus, 
we Search for the eyes only in the upper half of the face area, and the right eye and the left eye are searched within 
the right and left half of the region respectively. The mouth is searched only in the lower half of the face. After the 
eyes and the mouth are detected, the nose is searched within the region below the eyes and above the mouth. 

With the regions detected, the feature localization will be performed on the image patch of each specific region. For 
the eyes and the nose, the region image patches are cropped from the grayscale image, and for the mouth the region 
image patch is cropped from the color image. 


3 eye feature extraction 

3.1 overview 

We are interested in four feature points shown in Fig.2(a) as well as the eyelid contours. En (endocanthion) is the 
point at the inner commissure of the eye fissure, ex (ectocanthion) is the point at the outer commissure of the eye 
fissure, pi (palpebrale inferius) is the lowest point in the mid-portion of the free margin of the lower eyelid, and ps 
(palpebrale superius) is the highest point in the mid-portion of the free margin of the upper eyelid. These four points 
are usually of most interest to clinicians. In the following we describe the method for the right eye, and for the other 
eye we mirror the eye image first before the feature detection. 


3.2 iris detection 

The eye image is first smoothed by a Gaussian filter and normalized in intensity. The reflection in iris is then 
removed using the method in Asteriadis et al. (2009). Briefly, a binary image is generated using otsu’s histogram 
algorithm (Otsu 1979), and any white spot whose area is smaller than 20 pixels is filled with black (Fig.4(d)). After 
that, the corresponding pixels in the original image are filled with the intensity of their surrounding pixels (Fig.4(f)). 
Template matching method described in 2.1 is then used to locate the iris on the binary image, and the root template 
is a binary image of the iris. The detected iris region is shown in Fig.4(e). The top and bottom boundaries are 
extracted from the detected iris binary image. The top boundary will serve as part of the upper eyelid contour, which 
will be discussed in the next section. Two key points are located on the bottom boundary, which are the corners of 
the iris and the sclera (Fig.4(f)). They can be found by searching for the first curvature extrema from each side of the 
bottom boundary. 


3.3 eyelid extraction 

The eyelid model consists of an upper eyelid curve and a lower eyelid curve, connecting at corner points en and ex. 
As stated in Vezhnevets and Degtiareva (2003), edge map does not give a robust estimation of the eyelids, and 
luminance valley points along the horizontal lines were detected as the candidate upper eyelid points. Our method is 
mostly inspired by Vezhnevets and Degtiareva (2003). To detect the upper eyelid, a local minimum map is 
constructed where each local minimum point is the luminance valley point in either horizontal or vertical direction 
(Fig.5(a)). Since the middle part of the eyelid contour is already detected in 3.2, we only look for the upper eyelid 
points to the two sides of the iris. As shown in Fig.5(a), these two parts of the upper eyelid are well captured by the 
local minimum map, and two straight lines can approximately fit the local minimum points on these two parts. We 
apply Hough transform (Duda and Hart 1972) to detect the two straight lines from the local minimum map. Hough 
transform is known to be robust to imperfect data and noise. In order to reduce the search space, we restrict the 
locations of lines based on the iris detection result in 3.2. Particularly, the two key points in Fig.4(f) serve as the left 
and right bound for the right and left parts respectively, and a horizontal line passing the topmost point of the middle 
upper eyelid serves as the upper bound (Fig.5(b)). Hough transform produces a set of lines within each search range, 
and the line with maximum number of local minimum points is selected. The points on the selected lines are the 
candidate eyelid points (Fig.5(c)). 

Some false positive points that exceed the bound of the eye may be included in the detected lines. A corner 
localization procedure is applied to locate the eye corners from the candidate eyelid points. Only the eyelid points 
that are in between the two eye corners are retained, and others are discarded (Fig.5(d)). The template matching 
method in 2.1 is again applied for eye corner detection. The root templates here are the 9x9 patches centered at the 
left and right eye corners respectively. Instead of searching the entire eye image, we only go through the locations of 
each candidate eyelid point, and the point location at which there is a best match with the template is selected as the 
corner position. 

The retained points on the two sides of the upper eyelid together with the points on the middle eyelid curve detected 
in 3.2 constitute the upper eyelid points. A Bezier curve is fitted to all these points to obtain the upper eyelid contour 
(Fig.5(e)). Since the Bezier curve interpolates the end points, the resulting eyelid contour always pass through the 
two eye corner points. 

Fig.5(a) shows that the local minimum map does not capture the lower eyelid very well, because the contrast on the 
lower eyelid is not as strong. In Vezhnevets and Degtiareva (2003) and Zheng et al. (2005), a quadratic curve is 
fitted to the eye corners and the lowest point of the iris to generate the lower eyelid. However, this curve may not 
approximate the lower eyelid very well because a major portion of the lower eyelid in between the three points is 
missing. We approximate the lower eyelid by a piecewise linear curve which consists of four line segments. They 
are divided by the following five points: two eye corners (A,E), left lower bound of the iris (B), right lower bound of 
the iris (D), the lowest point in the middle of the iris (C) (Fig.6(a)). The two eye comers have been detected before 
so they are fixed, but the other three points are flexible to move. We start from the line segment connecting the left 
corner (A) and the left lower bound of the iris (B). By moving point B along vertical direction within a small range 
we can form several candidate lines. On each line the average gradient magnitude is calculated and the line with the 
maximum gradient magnitude is selected. The same method is used for the line segment between the right corner (E) 
and the right lower bound of the iris (D). Then we form several polyline candidates connecting the lower left (B), 
low middle (C), and lower right bound (D) of the iris, by moving the low middle point of the iris along the vertical 
direction. The polyline with maximum average gradient is selected. Fig.6(b) shows the points on the selected line 
segments for the lower eyelid. 


A Bezier curve is fitted to all the points on the above line segments to obtain the lower eyelid contour (Fig.6 (c)). 
The same as the upper eyelid contour, the lower eyelid contour also pass the two eye corner points. Therefore, the 
upper and lower eyelid contours connected by the same end points form an entire eye boundary (Fig.6(d)). 

With the upper and lower eyelid contour, The feature points in Fig.2(a) are easy to obtain. En and ex are already 
detected which are the two eye corners. pi and ps are the lowest and highest points on the eye boundary (Fig.6(d)). 


4 lip contour extraction 

4.1 overview 

Six feature points (Fig.2(b)) are of interest as well as the lip contour. Ch’s (cheilion) are the mouth corners, li 
(labiale inferius) is the middle point of the lower lip, Iss’s are the highest points on each side of the upper lip, and Is 
(labiale superius) is the lowest point in between the two Iss’s. These points are not only defined from clinical point 
of view, but also serve as parameters for the curve fitting which will be shown later. 


4.2 region based lip segmentation 

The mouth is rich in color so color images are often used for mouth feature detection. Various color spaces (HSV, 
Lab, Luv, RGB) have been explored (e.g.,Wang et al. 2004; Eveno et al. 2004; Yokogawa et al. 2007). A fuzzy C- 
means (FCM) method incorporating shape information (Wang et al. 2002) has been proposed to cluster the lip 
region based on selected color components. This method assumes an ellipse shape of the mouth and can better 
eliminate noisy regions produced from FCM. However, when the mouth is long and thin the ellipse model cannot 
approximate it very well. A new color transformation Eveno et al. (2001) represents each pixel as the curvature of a 
parabola defined by three points in the transformed color space. Generally this method can generate regions close to 
the mouth shape, but a threshold needs to be selected manually for different images. 

We apply a method similar to Lievin et al. (1999) which is much simpler than the above two methods. The ratio 
image R/G is computed, and then the pixel values of the resulting image are truncated to within 2.5 standard 
deviations around the mean pixel value. This can effectively remove the few extremely high ratios due to image 
noise. The histogram of this ratio image will have two peaks, the higher one being the skin and the lower one being 
the lip. Otsu’s method is applied again to find the optimal threshold on this ratio image. A binary image is then 
generated based on this threshold (Fig.7(b)). Connected components algorithm (Dillencourt et al. 1992) is performed 
to remove small background noise and holes in the mouth region. Then the boundary of the mouth is extracted to get 
an initial contour (Fig.7(c)). 


4.3 parabola model for the lip contour 
We use parabolas to model the lip contour. Four parabolic curves (yl, ylr, yrl, yr) are used to approximate the upper 


lip and one (yb) to approximate the lower lip (Fig.8). The intersection points of the upper lip curves are Dis D, and 
D3, where Pp, and Pp, are also the zero-derivative points of the connecting parabolas. P, is the zero-derivative 
point of curve yb. The reason to use four parabolas instead of two for the upper lip is that the curve between D, and 
DP may not follow the same parametric form as the curve on the other side of Pp, and p,. Moreover, bending and 
expansion are applied to the point (xX, y) on the curve as follows. 

x'= T(B(x, y), y) (5) 

where X' is the new x coordinate of point (x, y) and the y coordinate does not change. B(x, y) is the bending 
function (Wo” rz and Rohr 2006) and T(x, y) is the expansion function defined as follows: 


B(x, y)=x—y?d , T(x, y) =(x-x,)/exp(y-s)+X (6) 
where X, is the x coordinate of the zero-derivative point of the parabola and 0, S are coefficients. The sign of s 


determines the shrinking or expansion of the curve along x direction. 
The origin of the coordinate system is at the left mouth corner, and the x,y axes are shown in Fig.8. The coefficients 


of each parabola curve can be derived from points P,~ P, as well as the two mouth corner points Pp, and p,. The 
parameters to be optimized are the coordinates of these six control points and the coefficients for expansion and 
bending of each curve. 


The optimization is performed by maximizing an energy function. The energy function in (Wang et al. 2004) is 
based on pure region information obtained from clustering. If there are some misclassified pixels in the probability 


map, the fitted curve will likely copy this error. We use an energy function as follows which combines region and 
gradient information. 
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where grad (p) is the normalized gradient magnitude at point p , h(k,I) is the kth histogram bin value of image I, 


(7) 


Q,, and oom are the regions inside and outside the mouth in a small neighborhood of point p , and C includes the 


entire lip contour. region( Pp) defines the region difference in the neighborhood of each point. @ is the coefficient 
to balance the two terms. 


4.4 parameter initialization and optimization 

The rough positions of the control points can be estimated from the initial mouth contour obtained in 4.2. First, the 
leftmost and rightmost points are identified as mouth comers. After the two corner points are located, we detect the 
rough locations of other four control points from the initial contour. The coordinate system in Fig.8 is first 
constructed, with the line connecting the two corner points as the x axis. Under this coordinate system, the middle 
point in the x direction on the upper lip is first located, and the point in the vicinity of this point with the highest y 


coordinate is Ds (note that y axis is pointing down). D, divides the upper lip into left and right parts. On each part 
the point with smallest y coordinate is Di and Ds respectively. On the lower lip contour the middle point in the x 


direction is p 4- Fig.9(a) shows the estimated initial positions of the control points. 
Other parameters to be optimized are the bending parameter 6 and the expansion parameter S for each parabola 
curve. We set S=0 and 6 =0.0001 for each curve. The optimization is constrained by the following 


inequalities to ensure reasonable spatial relation among P,(X,, Y,), D>(X2,Y¥2) and P3(X3,Y3)- 


My SAS Mao? Vie 7 V3 (8) 

Moreover, the mouth corners are constrained on the minimum luminance line, so only the x coordinates of the moth 
corners need to be optimized. 

The feature points in Fig.2(b) coincide with the control points in Fig.8, so the optimized parameters of the parabola 
curves give the coordinates of the feature points as well. Fig.9(b) shows the lip contour and the feature points after 
optimization. 


5 nose boundary extraction 

5.1 overview 

Very few literatures have reported feature detection on the nose, probably because the nose does not have distinct 
feature points as the mouth and the eyes do. In this paper, we are interested in the alare point which is the most 
lateral point on each side of the nose (Fig.2(c)), and the nose boundary passing through the alare point on each side. 
However, since each boundary is an open contour, there is no clear upper or lower bound of the contour. 


5.2 nostril detection 

In order to facilitate the detection of the nose boundary, we first locate the rough positions of the nostrils. The 
average intensities of the rows and the columns of the image are calculated, generating a horizontal profile and a 
vertical profile of the intensity (Fig.10(a)). The x coordinates of the nostrils can be found at the local minima on the 
horizontal profile, and the y coordinate can be found at the local minimum on the vertical profile The search range 
of the nose boundaries is then limited to the two sides of the nostril points (Fig.10(b)). 


5.3 nose boundary detection 

Since the nose boundary is usually weak, histogram equalization is first applied on the grayscale nose image to 
enhance the contrast (Fig.11(a)). An edge map of this enhanced nose image can usually capture most of the nose 
boundaries. However, it is not always possible to get perfect edges around the nose boundaries, and an example of 
the missing edge around one boundary is shown in Fig.11(b). Alternatively, local minimum map (Fig.11(c)) can be 
applied here as well which is very close to the edge map around the boundary points, but similar problems may 
occur as with the edge map. Therefore, we decide to fuse the edge map and local minimum map into a combined 
map that can capture the boundary points more robustly. Specifically, a local minimum point is added to the edge 


map only if there is no surrounding edge point in a small neighborhood of this local minimum point. Then the 
edge/local minimum points (or pixels) are grouped by connected components algorithm, and small components are 
discarded. Fig.11(d) shows the combined map from Fig.11(b) and (c) after small components removal, where the 
missing edge is filled in by local minimum points. 

Generalized Hough transform for arbitrary shapes is then used to detect the boundary contours from the combined 
map. A template shape and a reference point need to be defined first, and a lookup table is constructed recording the 
location of each point on the template shape relative to the reference point. Given the location of the reference point, 
a shape can be generated from the lookup table, thus each point in the Hough space (ie., reference point) 
corresponds to a shape in the image space. Assuming that the rotation and scaling of the shape do not change, we 
have the coordinates of the reference point as the only parameters. Since these parameters are equivalent to those in 
the line detection with Hough transform, the algorithm for the curve detection follows the same routine as the line 
detection. In our implementation, the template for each nose boundary is a generic shape of the nose boundary, and 
the reference point is chosen as the leftmost (rightmost) point for the left (right) boundary (Fig.11(e)). With the 
results of the nostril locations, the range of the parameters is restricted to the part of the image shown in Fig.10(b). 
Rotation and scaling parameters can also be included, but the computational cost will be significantly increased. We 
find that a fixed template shape can serve our purpose as well since the nose boundaries do not vary too much 
among different individuals, so we do not consider rotation and scaling parameters. 

A Bezier curve is then used to fit the detected local minimum points on each boundary. The alare points are then 
located by finding the left most and rightmost point on each contour respectively (Fig.11(f)). 


6 results and validation 

The proposed method is tested on FERET database (Phillips et al. 1998; Phillips et al. 2000), which contains various 
color images of faces with different expressions, rotations and lighting conditions. From the frontal face subset, the 
images that obviously violate our assumptions are removed, and 100 images are randomly selected from the 
remaining images which meet our assumptions to a maximum extent. Slight facial expressions and head rotations 
are tolerable for our method. Each image has the dimension of 512*768 pixels and includes the full head and the 
neck of a single person against a near-white background (Fig.3(a)). Fig.12 shows several results on the detected face 
images, from which we can see the extracted feature points and boundary contours are faithful to the image data. 

In order to quantitatively evaluate the performance of our method, the feature points are manually labeled on each 
image by three different operators. The ground truth position (P) is taken as the average of the manually labeled 
positions from the three operators. For each feature point on an individual image, the accuracy of our detection 
result is calculated as 


a=1-||P,-P,|l/4, © 
where D gis the ground truth position of this feature point, p, is the position obtained by our method, ||.|| is the 


Euclidean distance, and d_, is the width of the corresponding facial component. For example, d_, is the eye width 


when the accuracy of an eye feature point is calculated, which can be estimated from the distance between the 
ground truth positions of two eye corners. The value of a is truncated to [0,1]. The result of a feature point is 
considered successful if its accuracy is above a threshold T, and the result of an entire image is successful if all the 
feature points are successfully detected. The overall success rate for each feature point is the percentage of the 
successful cases, and the overall success rate for all feature points is the percentage of successful images. 

Among all the results, there is no failure in face detection or region detection. Note that the face detection is 
considered successful if all the facial components of interest are completely inside the detected face, and the region 
detection is considered successful if all the related feature points are inside the detected region. Among all the 
failure cases, the main cause of the mouth failure is the incorrect separation between the lip and the skin due to the 
lack of color constancy under poor lighting. The failure on the eye features is mostly because the iris is connected to 
other parts in the binary image thus leading to incorrect iris detection. The failure of the nose features is mainly due 
to extremely low contrast on that particular boundary. Table 1 shows the average accuracy of each feature point, the 
success rate for each feature point and the overall success rate. The success rate of individual feature point remains 
high even with a high threshold T, although the overall success rate declines faster with the increasing T because it 
is affected by every single feature point. 

No validation of the boundary contours is performed because it is difficult to provide the ground truth by manually 
delineating the contour. Since most of the feature points are obtained from the contours, the accuracy of the feature 
point detection can reflect the accuracy of the boundary extraction to some extent. 


7 conclusion 

This paper proposes a framework for automatic facial feature detection. Features of the eyes, the nose and the mouth 
are extracted. To our knowledge this is the first work to extract a complete set of major facial feature points and 
boundary contours. The primary application of this work is the computer aided diagnosis, so only frontal images are 
considered and some other assumptions of the input image are made. Visual inspection of the experimental results 
shows the accuracy of the overall feature detection, and quantitative validation further demonstrates the accuracy of 
the feature point detection. 
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Figure legends: 

Fig.1 framework of the facial feature detection 

Fig.2 feature points and boundary contours of the (a) eye (b) mouth (c) nose 

Fig.3 (a) the input color image (b) the image converted to grayscale and cropped to the detected face (c) template of 
each facial component 

Fig.4 (a) the original eye image (b) the result of Gaussian smoothing and intensity normalization of (a) (c) the 
binary image of (b) (d) the binary image after hole filling (e) the detected iris region (f) boundaries and key points 
extracted from the iris region, shown on the eye image after reflection removal 

Fig.5 (a) local minimum points overlaid on the eye image (b) search range of the upper eyelid points (the lower left 
and lower right rectangle regions) (c) detected upper eyelid points by Hough transform (d) the entire set of upper 
eyelid points after eye corner detection (e) final upper eyelid contour 

Fig.6 (a) piecewise linear curve of the lower eyelid (b) detected lower eyelid points (c) final lower eyelid contour (d) 
the entire eyelid boundary and the feature points 

Fig.7 (a) the original mouth image (b) the binary image (c) the mouth contour extracted from (b) 

Fig.8. pl, pr: mouth corner points; p1~p4: four points used to determine the parameters of the parabolas; y*: 
parabolas for different segments of the lip 

Fig.9 (a) initial lip contour and the minimum luminance line (green line) (b) initial control points (c) final lip 
boundary and feature points 

Fig.10 (a) illustration of the nostril detection (b) detected nostrils and the search range of the nose boundaries 
(rectangle regions on the left and right sides) 


Fig.11 (a) nose image after histogram equalization (b) edge map of (a) (c) local minimum map of (a) (d) combined 
map after small components removal (e) template nose boundaries and the reference points (f) detected nose 
boundaries and feature points shown on the original nose image 

Fig.12 detected feature points and contours on several color images 


Figures: 
Eye feature extraction 
Face detection Region detection 
mouth feature extraction 
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Table 1 success rate of each feature point and the overall success rate 


Facial components | Feature points | Accuracy Success rate 
T=0.95 T=0.9 

Left eye ex 0.952 0.68 1 

en 0.943 0.63 1 

ps 0.919 0.61 0.75 

pi 0.945 0.72 0.89 
Right eye en 0.941 0.62 1 

ex 0.938 0.68 1 

ps 0.912 0.68 0.77 

pi 0.948 0.70 0.91 
Nose Al (left) 0.903 0.58 0.63 

Al (right) 0.938 0.75 0.89 
mouth Ch (left) 0.993 1 1 

Lss (left) 0.988 0.99 0.99 

Ls 0.988 0.99 0.99 

Lss (right) 0.984 0.99 0.99 

Ch (right) 0.963 0.84 0.95 

Li 0.959 0.88 0.96 
Overall 0.52 0.90 


,D Shape Analysis of Lateral Ventricles in Autism 


= Duan *», Judith Miles” 


der of unknown etiology 
ments. The ventricular system 


in autism has been reported in previous works. Powever, the magnitude of the 
differences revealed by volume analysis is limited, because it does not consider how 
shapes overlap each other. Shape analysis can precisely locate morphological changes 
in brain structures which cannot be reflected in volume measurements, thus it 
becomes more and more popular in neuroimaging community. Very few studies have 
been done with regard to the shape morphology of the ventricles in autism. 


The aim of this paper is to study the shape differences of the lateral ventricles 
between autistic children and normal controls. 3D models of ventricle shapes are 
compared to reveal the shape difference at every surface location. 


Participants 


Five Children with autism were recruited from the Thompson Center for Autism and 
Neurodevelopmental Disorders. Four control subjects with matching gender, age and 
ancestry were recruited from the community. Student t-test was used to compare the 
ages and x2 test was used to compare the gender ratios. There was no significant 
difference between the two groups in terms of age, gender and race. This study was 
approved by the Health Sciences Institutional Review Board. The parents or legal 
guardians of all subjects provided written consent for participation in this study. 


Methods 


The 3D ventricle models were obtained using the semiautomatic segmentation 
software that has been developed in our previous work. Six anatomical landmarks 
were identified. Pseudo-landmarks on the surface were then interpolated using the 
method in our previous work. Generalized Procrustes Analysis is performed on all the 
3D shapes in order to eliminate the shape differences caused by translation, rotation 
and scaling. difference between groups at every surface location was tested using 
Hotelling T2 two-sample metric. Since comparisons are made at thousands of surface 
points, we adopted False Discovery Rate estimation (FDR) for p-value correction. 
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Results 

From the overlaid average structures of the two groups (Fig.1(a)), we find that the 
posterior horns of both ventricles in autistic children are bent more inward, and the 
anterior horns are bent more downward in the autistic children. The amount of 
difference can be seen from the signed distance map (Fig.1(b)). A two-tailed alpha level 
of 0.05 is chosen as the significance threshold for the raw p-values. The raw significance 
maps of both ventricles show significant shape difference in part of the body and the 
posterior horns between autistic children and controls (Fig.1(c)), but the corrected 
significance maps show no shape difference between the two _rou, s (not shown). 


Conclusion 

This paper studies the shape differences of the lateral ventricles between autistic 
children and controls. The overlaid mean structures show shape difference in the 
anterior and posterior horns, but sstatistical tests do not show any significance in those 
differences. This may be due to small sample size and further experiment needs to be 
done when more data are available. 


Left ventricle Right ventricle 


p>=0.05 (b) ] 
™. 


Fig.1 on each side: (a) the average shape of the 
pathological (red) and normal structure (blue) (b) 
signed distance map (negative distance indicates the 
pathological shape is inside) (c) raw significance map 
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SUMMARY 


A number of studies have documented that autism has a neurobiological basis, but the anatomical extent 
of these neurobiological abnormalities is largely unknown. In this paper, we apply advanced computational 
techniques to extract 3D models of the corpus callosum (CC) and subsequently analyze local shape 
variations in a homogeneous group of autistic children. Besides the traditional volumetric analysis, we 
explore additional phenotypic traits based on the oriented bounding rectangle of the CC. In shape analysis, 
a new conformal parameterization is applied in our shape analysis work, which maps the surface onto a 
planar domain. Surface matching among different individual meshes is achieved by aligning the planar 
domains of individual meshes. Shape differences of the CC between autistic patients and the controls are 
computed using Hotelling T? two-sample metric followed by a permutation test. The raw and corrected 
p-values are shown in the results. Additional visualization of the group difference is provided via mean 
difference map. Copyright © 2009 John Wiley & Sons, Ltd. 


Received 7 February 2008; Revised 24 April 2008; Accepted 19 February 2009 


KEY WORDS: corpus callosum; autism; shape analysis; planar conformal mapping; permutation test; 
phenotypic traits 


1. INTRODUCTION 


A number of studies have documented that autism has a neurobiological basis, but the anatomical 
extent of these neurobiological abnormalities is largely unknown [1]. Several studies have reported 
deficits in the size of the corpus callosum (CC) and its sub-regions, although the results are 
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inconsistent with regard to CC sub-regions. For example, Piven et al. [2] reported reductions 
in the size of the body and posterior regions of the CC in autistic patients, Hardan et al. [3] 
found significant differences in anterior regions, and Vidal et al. [1] found reductions in genu 
and splenium as well as the total CC area. A detailed summary of CC abnormalities in autism 
can be found in [4]. The inconsistency of the results may be due to factors such as the sample 
size, subject age and gender, and specific diagnosis. Particularly, heterogeneity within the autism 
diagnosis obscures the genetic basis of the disorder [5]. Recently, Miles et al. proposed a new 
definition of autism subgroups, which divided autism into essential autism and complex autism [5]. 
Limiting studies of brain morphology to individuals with essential autism decrease the background 
noise of structural variation and allow the analysis of the more uniform population [5]. He ef al. 
[6] first adopted this definition and compared the shape difference between subjects with essential 
autism and controls, although their sample size was conspicuously small. In this paper, we also 
focus on essential autism, which will greatly reduce the heterogeneous factors. 

Quantitative morphologic assessment of individual brain structures is often based on volumetric 
measurements and shape analysis. Shape analysis has gained an increasing amount of interest 
from the neuroimaging community because it can precisely locate morphological changes in the 
pathological structures, which cannot be reflected in volume measurements. A significance map 
is often computed in shape analysis, which tells how significant the difference is between two 
groups of structures at each location. However, shape analysis gives little description with respect 
to the phenotype, which is directly linked to the genetic basis of the diseases. A phenotype is the 
appearance of visible characteristics of an organism produced by its genes and their interaction 
with the environment; a phenotypic trait is a visual category of phenotypic variation [7]. Under 
this definition, CC volume is a trait, but agenesis of the CC is a phenotype. Most of the previous 
works only focus on the volume measurement. However, there are many other traits that can be 
defined based on the characteristics of the CC shape. 

Shape morphology of the CC in autism has been studied in [1,6]. In [1], CC thickness (the 
distance from a medial line) at each surface point was compared between patients and controls. 
In [6], an average CC model was used as a template and the distance to the template at each point 
was compared between two groups. However, the shape analysis in both was still in a 2D manner, 
because the comparison was performed on a contour-by-contour basis. 

This paper gives a comprehensive pipeline of the analysis of the CC with respect to the shape, 
morphology and phenotypic traits in essential autism (Figure 1). The 3D model of the CC is 
reconstructed from 2D contours of nine sagittal slices as in [1]. We use a newly developed 
semiautomatic method [8] to segment the CC from 2D MR images. Compared with the manual 
tracing method in most previous work, our method is faster and more accurate, which will facilitate 
the following analysis. In the phenotypic trait analysis, volumes of the CC and five sub-regions are 
compared in a way similar to [1]. CC thickness at each region border is also compared. Moreover, 
an oriented bounding rectangle of the 2D CC on sagittal magnetic resonance imaging (MRI) is 
constructed, and the length, width and ratio (length/width) of the bounding rectangle are the MRI 
traits to be compared between autism and control groups. Since the oriented bounding rectangle 
of an object is the minimum rectangle that encloses this object, the measurements of the bounding 
rectangle can give both size and shape information. In the shape analysis, we directly compare the 
3D coordinates at each surface point of the CC model using Hotelling T* two sample tests [9]. In 
order to perform this comparison, point correspondence among different individual models must be 
established. We develop a new planar conformal mapping based on the methods in [10, 11] to map 
the surface onto a planar domain and thus find a one-to-one point correspondence among different 
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Figure 1. Schematic view of the entire procedure. 


models. The results of this paper show both the differences in the phenotypic traits measurements 
and the visualization of local shape differences between patients and controls. 


2. METHODS 


This section describes the details of our methods, including CC segmentation, phenotypic trait 
measurement and statistical shape analysis. The entire procedure can be summarized in the 
following steps: 


1. For each subject, CC contours are segmented from mid-sagittal slice and four adjacent slices 
on each side, and a 3D model of the CC is constructed by contour stitching. 

2. CC volumes and measurements of the oriented bounding box are computed, and differences 
of these measurements are compared between patients and controls using f-tests. 

3. Each CC model is mapped onto a planar domain and re-triangulated to establish the point 
correspondence among different models. 

4. Spatial alignment of all CC models is performed using Procrustes analysis. 

5. 3D coordinates at each surface point are compared between patients and controls using 
Hotelling T2 tests, which results in a significance map of group differences. 


The following subsections will describe each step in details. 


2.1. CC modeling 


We start with slice-by-slice segmentation and stack the 2D curves to make a 3D model. Because, 
it is more straightforward to verify the accuracy of 2D results slice-by-slice, this method provides 
better results as opposed to direct 3D segmentation and validation. Mid-sagittal slice and four 
adjacent slices on each side are used for CC segmentation. The local shape variation of the CC is 
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fornix 


splenium 
rostrum 


Figure 2. Four parts of the CC boundary: anterior (AC), upper body (CD), posterior (DB), lower body 
(BA) (image modified based on original illustration from [12]). 


(d) 


{e) a 


Figure 3. (a) User-initialized polyline; (b) the edge map; (c) the completed seed contour; (d) distinct 
points around the fornix; (e) the result without fornix removal; and (f) the result with fornix removal. 


dramatic, and parts of the CC may be narrow or have bumps. The most challenging problem is the 
existence of the fornix, which is a thin structure that may or may not contact the CC in sagittal MR 
images; it is almost the same intensity as the CC. He et al. proposed a part-based semiautomatic 
method [8] to segment the CC on sagittal slices, which can effectively overcome most of the 
challenges. The novelty of this method is dividing the CC boundary into four parts connected by 
four sensor points (Figure 2). The user initializes a polyline inside the CC region by three mouse 
clicks (Figure 3(a)), and a seed contour consisting of four parts is automatically generated by point 
tracing on the edge map (Figure 3(b) and (c)). The details of the seed construction can be found 
in [8]. The seed will then automatically deform according to active contour evolution, but each 
part has its own motion law. The contour evolution mechanism can be found in [8]. 

The fornix tip can also be removed by automatic fornix detection [8]. Since we know the fornix 
always appears beneath the body of the CC, we only need to search along the lower boundary to 
detect it. There are three distinct points around the fornix (Figure 3(d)). The fornix dip is removed 
by connecting points a and b. The final result on one sagittal slice is shown in Figure 3(f), compared 
with the intermediate result with the fornix dip in Figure 3(e). The quantitative evaluation of the 
segmentation results is performed in [8], which shows high accuracy of this method. 

Since the CC shapes on adjacent slices do not differ too much, we apply the segmentation 
method on mid-sagittal slice and the result is used as the seed for its two neighbor slices. After 
slight deformation, the boundary curves on these two slices can be obtained and each one is again 
used as the seed of its next adjacent slice. 

Contour stitching is performed to create the 3D CC model of each subject. Since the segmentation 
method can keep track of the four sensor points on the CC boundary, the point correspondence 
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Figure 4. (a) Partition of the CC: rl—anterior third; r2—anterior body; r3—posterior body; r4—isthmus; 
r5—splenium and (b) bounding rectangle of the CC. 


problem can be simplified by first matching the sensor points among contours [6]. After matching 
the sensor points, the four parts are implicitly matched among contours. We only need to match 
the individual points within each part. The details of matching the sensor points and the inner 
points of the four parts can be found in [6]. 

By connecting the matched point pairs between adjacent contours, the 3D mesh is created. 
A tangential Laplacian smoothing [13] is performed to maintain a good node distribution of the 
model. The reconstructed model is shown in Figure 5. 


2.2. Phenotypic traits 


Volumes of the CC and its sub-regions are computed similar to that in [1]. The CC is divided 
into five regions along anterior—posterior line (Figure 4(a)). To find the anterior—posterior line, we 
calculate the distance between every two points on the CC boundary, and find the pair of points 
with the longest distance (pl and p2 in Figure 4(a)). The line connecting these two points is the 
anterior—posterior line. The five regions are defined the same as in [1]. For simplicity, we label 
them from rl to r5 as shown in Figure 4(a). Areas are computed in pixels for each 2D segmented 
CC, and the areas of nine slices are summed to generate the voxel count of the 3D model. The raw 
volume is the multiplication of the voxel count and the voxel size in mm?. To take into account 
the effect of the brain volume, we normalize the raw CC volume by the total brain volume (TBV) 
and the intracranial volume (ICV), respectively. TBV includes gray matter and white matter and 
excludes ventricles [14], and ICV is the sum of white matter, gray matter, and inner and outer 
cerebrospinal fluid spaces [15]. A choice can be made between using TBV or ICV as an adjustment 
factor [14], but our results show that the two choices do not make any difference. The brain 
volumes, raw and scaled volumes of the CC and the five regions are compared between patients 
and controls using t-tests. 
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We explore additional phenotypic traits based on the oriented bounding rectangle of the CC. To 
construct the rectangle, we first find points pl and p2 as mentioned above. On each side of line 
plp2, we find a point on the CC boundary which has the maximum distance to line pl p2 (p3, p4 
in Figure 4(b)). The bounding box is a rectangle whose long edges pass p3 and p4, respectively 
and are parallel to line plp2, and whose short edges pass p1 and p2, respectively (Figure 4(b)). We 
measure the length, width and aspect ratio (length/width) of the bounding rectangle on each 2D 
slice and average them across the nine slices. The shape of the bounding box depends on both the 
size and the shape of the CC, so these traits can reflect both size and shape information of the CC. 
We also compute the thickness at the border of every two sub-regions, which is the length of each 
dividing line (L1-L4 in Figure 4(a)). With a little abuse of notation, we denote the thickness at 
each dividing border as L1—L4. The thickness at each border is averaged across the nine slices. The 
average measurements of the bounding box and the thickness are compared between patients and 
controls using t-tests. To account for the brain volume effect, we scale the above raw measurements 
by the cubic root of TBV and ICV, respectively. Since these measurements are in millimeter, this 
scaling makes the units consistent. These scaled measurements are also compared between patients 
and controls. 


2.3. Point correspondence via planar conformal mapping 


Conformal parameterization has been explored intensively as a potential approach to the matching 
and analysis of brain data. It maps the brain surface into regular and more simple domains and 
carries out the analysis on the parameter domain. Gu ef al. [11] used spherical conformal map 
to match cortical surfaces, where they parameterized the genus zero closed surfaces by spheres. 
Unlike their data, the CC surfaces we are studying have two boundaries, which are the 2D CC 
contours of the two end slices. Inspired by [10], we design a planar conformal mapping method. 
In order to find point correspondence among different CC surfaces, we first parameterize each 
CC surface using a planar domain, and then align the parameterization in the planar domain so 
that surface points across different objects posses the same parameterization. Finally, the aligned 
parameterization in the planar domain is mapped back onto the original surfaces (Figure 5). 
Among all the CC models an arbitrary model is chosen as the target surface, and the corre- 
sponding points on the other models (source surfaces) are to be established to match the points on 
the target surface. We conformally flatten the target surface M onto a planar domain M’, which 
is a rectangle with the two boundaries y, and y, mapped to the upper and lower long edges. 
In order to parameterize the given surface with a planar domain, we need to construct a special 


Figure 5. Illustration of establishing point correspondence via planar conformal mapping: two individual 
CC models (left) and their planar domains (right). 
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holomorphic one-form t+7./—I, such that its integration around the boundaries has a real part 
with a constant period, while the integration along a path from one boundary to the other has a 
constant imaginary part. To get the imaginary part 7, we need to compute a harmonic function f 
by solving a Dirichlet problem: 


Af =0 
fly, =1 (1) 
fly, =0 


n is the gradient of f and the real part t is the conjugate of y. In our implementation, we first 
convert 7 into a piece-wise constant vector field v, sampled at each face; rotate v, by 90° around 
the face normal, and turn the rotated vector field into the one-form t. Choosing an appropriate 
starting point, we can integrate the holomorphic one-form t+7/—1 on M. Taking the real and 
imaginary parts of the integration to be y and x coordinates, respectively, we get a horizontal strip 
that is periodic along the X direction. Cutting the strip along the vertical line x =0, we get the 
final planar domain M’, where the upper and lower side corresponds to the two boundaries in the 
original mesh M. The map from M to M’ is conformal; the height—width ratio of the rectangle is 
a conformal invariant of the original surface. In our experiments, we always re-scale the rectangle 
uniformly along the X and Y directions so that the width equals to 27. As a consequence, the 
intrinsic structure is solely reflected by the height of the rectangle for our cases. 

The same method is used to map each source surface onto its planar domain. After that, we 
compute the map between the target planar domain and each source planar domain. In our case, a 
reasonable map should satisfy two constraints. First, each boundary of the source surface should 
be aligned with the corresponding boundary of the target surface. In the rectangle domain, we 
need to map the upper (lower) edge of the source rectangle to the upper (lower) side of the target 
rectangle. As mentioned before, the width of the rectangle is always equal to 27, while the height 
may be different from model to model. Therefore, we may need to re-scale one rectangle along 
the Y direction to align with the other rectangle. In our experiments, all the rectangles have similar 
heights, so the re-scaling factor is very close to one and therefore would not hurt the conformality 
of the map too much. Even in the case of a bigger re-scaling factor, the map is still affine and 
therefore harmonic, which is already strong enough for our purpose. The second constraint is 
to align the ‘bending corner’ (i.e. point c in Figure 2(d)) across two surfaces. In order to make 
the mapping process easier, we choose the bending corner to be the starting point used in the 
integration of the holomorphic one-form. As a consequence, in the planar domain the bending 
corner will appear at the vertical cutting side. Aligning the two rectangles along their left and right 
sides, the bending corners will be automatically matched. Since all the rectangles have the same 
width 27, we do not need to re-scale any of them along the X direction. 

Once the source rectangle is totally aligned with the target one, we get a map F’ between two 
planar domains. For each vertex v on the target rectangle, there is a point F’(v) in the source 
rectangle, which may lie inside a triangle, on an edge or at a vertex. Without loss of generality, 
we represent the point as a linear combination F’(v) =bov9 +b, v1 +b2v2, where v; is the vertex 
of the resident face and b; the corresponding barycentric coordinate. The map F’ between two 
planar domains induces a map F between the two original surfaces. Namely, each vertex v in the 
target surface corresponds to a point F(v) =bov9 +b1v1 +202 in the source surface. To facilitate 
later processing, we re-sample the source surface using these points as vertices, and re-triangulate 
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the source surface with the same connectivity as that of the target surface. Therefore, each point 
on the source surface has a correspondence on the target surface, and the two surfaces have the 
same triangulation. 

The same method is used to align each source planar domain with the target planar domain, and 
to re-triangulate the source surface. In this way, the point correspondence among all the surfaces 
is established through the target surface. 


2.4. Spatial alignment 


Before the shape comparison, rigid-body Procrustes Analysis [16] is performed on each 3D mesh 
across subjects in order to eliminate the shape differences caused by translation, rotation and 
scaling. The algorithm can be summarized as follows: 


1. Select one shape to be the approximate mean shape (e.g. the shape of the first subject). 
2. Align each shape with the approximate mean shape: 


(a) Translate each shape so that its centroid is at the origin. 

(b) Normalize each shape X by a scale factor s=||X||, where X is the coordinate matrix 
consisting of row vectors of the coordinate of each point after translation. 

(c) Rotate each shape to align with the newest approximate mean shape. The details of the 
rotation algorithm can be found in [16]. 


3. Calculate the new approximate mean shape from the aligned shapes. If it is different from 
the one in step 2, return to step 2. 


2.5. Statistical testing of group mean differences 


The difference in the size of the CC has been eliminated in the spatial alignment, so the shape 
analysis only reveals pure shape difference in local CC area between patients and controls. In the 
local shape analysis, difference between groups at every surface location is tested (Figure 6). This 
can be done in two main ways [9]. One way is to analyze the local difference to a template, which 
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Figure 6. Schematic view of the test procedure. 
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is usually the mean of the two groups. Student f-test can be used to test the local significance and 
He et al. [6] used this method. The main disadvantage of this method is the need to select a 
template, which introduces an additional bias. The other method is to directly analyze the spatial 
location of each point, and we apply this method in our study. No template is needed in this 
approach and multivariate statistics of the (x,y,z) location are used. Hotelling 7? two-sample 
metric is used to measure how two groups are locally different from each other. We use a modified 
T* metric instead of the standard one, because it is less sensitive to group differences of the 
covariance matrixes and the sample size [9]. This metric is defined as 


-1 
T? = (14 —py) (ca/my5+C/na) a) (Hy — My) (2) 


where }°, and >", are the covariance matrices of the two groups. Since comparisons are made at 
thousands of CC surface points, a permutation test is, therefore, used to confirm the significance 
of the overall differences in the statistical mapping result adjusting for multiple comparisons [17]. 
Suppose we have two groups A and B whose sample sizes are na and ng, the permutation scheme 
can be stated as follows: 


1. Calculate the test statistic (in our case, T?) between A and B, denoted as So. 

2. The observations of A and B are pooled. From the pooled values, 7 4 observations are sampled 
without replacement to form group A, and ng observations are sampled to form group B. 
The test statistic between the two new groups is calculated. 

3. Repeat the sampling process (permutation) in step 2 for M times. The test statistic in jth 
permutation is denoted as $j. 

4. The proportion of S;’s, which are greater than So is the p-value. 


M is a large number with the maximum value of 


natnp 
NB 


The resulting p-values generate a significance map that locates significant shape differences 
between the two groups. Both raw and corrected p-values are shown in the experiment. 

We also perform a global shape analysis, where the average of the local Hotelling T* metrics 
across the whole surface is calculated for each subject, and the p-value is computed in the same 
way as for the local testing. No permutation is needed since there is only one test. 


3. RESULTS 


3.1. Subjects and data 


Thirty patients with essential autism participated in this study. All patients met the DSM-IV 
criteria [18] as well as autism diagnostic interview-revised (ADI-R) and ADOS algorithm criteria. 
The distinction between essential and complex autism is based on the presence of generalized 
dysmorphology, which is indicative of an insult to early embryologic development and has been 
previously reported [5]. Children with a known disorder such as Fragile X, Tuberous Sclerosis, 
a chromosome abnormality or severe prematurity with brain damage and children with IQs/DQs 
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Table I. Volume measures. 


CC measure Patients Controls p-value 
TBV (cm?) 1421.94+299.4 1313.6+336.0 0.09 
ICV (cm) 1460.44 308.3 1348.3+345.5 0.09 
Raw measures (mm?) 

Total CC 4229.1+4755.3 4727.7.04971.4 0.01 
rl 1754.2+321.3 1939.5 +367.9 0.02 
12 474.4+ 103.9 556.14 142.9 0.008 
13 439.74 123.1 505.14139.9 0.03 
14 404.6 + 108.3 458.8 + 140.8 0.05 
15 1153.34227.7 1263.5+275.7 0.05 
Scaled by TBV (x 1074) 

Total CC 31.049.3 37.0+8.1 0.004 
rl 13.043.5 15.143.1 0.004 
12 3.4+1.0 4.3+1.0 0.001 
13 3.21.2 3.9+1.0 0.009 
14 2.9+1.0 3.6+0.9 0.01 


si) 8.642.9 10.142.4 0.02 


less than 40 will be excluded. Five healthy people recruited from the NIH Human Brain Project 
participated as controls. Besides, we also obtained MR images of 23 controls from IBSR database 
[19]. The age range of all participants is 3-30. 

MR images of our patients and controls were obtained from our autism aesthesia protocol, 
which included intravenous infusion of protocol to insure image quality. Axial, coronal and sagittal 
Tl-weighted images were acquired using the Siemens Symphony 1.5 T scanner with the following 
parameters: TR=35ms, TE=minimumms, NEX=1, flip-angle=30°, thickness=1.5mm, field 
of view =22cm, matrix =512 x 512. The data obtained from IBSR were T1-weighted MR images 
with 3.1 mm thickness. Each brain volume was interpolated to 0.9mm isotropic voxels. 


3.2. Phenotypic trait measurements 


Table I shows the t-test results of the volumes. The mean TBV and ICV of the patients are 
greater than those of the controls, but the p-values do not reveal any significance in this difference 
(p>0.05). The raw measures of the CC volume and its sub-regions show that the patient’s volume 
is significantly smaller than the controls (p<=0.05), especially in the region r2 (anterior body). 
The comparison of the scaled measures by TBV and by ICV shows similar results, both of which 
augment the significance of the difference. Table I only shows the results of TBV normalization. 

Table II gives the t-test results of the traits regarding the bounding rectangle and the thickness. 
In the raw measures, the patients have significantly reduced length and aspect ratio of the bounding 
rectangle, while other traits do not show any significance in the difference between patients and 
controls. When scaled by the cubic root of TBV, the length of the bounding rectangle keeps its 
strong significance and the width still shows no significance in the group difference. The p-values 
of the difference in thickness are decreased but none of them reach the significance level (0.05). 
The difference in L2 between patients and controls is close to significant, which is consistent 
with the volume measurements since the anterior body displays most significant difference. The 
comparison shows similar results when the measurements are scaled by the cubic root of ICV, and 
Table II only shows the results of TBV normalization. 
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Table II. Phenotype traits. 


CC measure Patients Controls p-value 
Length 66.2+4.8 72.145.0 <10-4 
Width 24.0+3.0 23.6+3.6 0.33 
Length/width 2.78+0.3 3.11+40.49 0.002 
Ll 5.914+1.13 6.09+ 1.19 0.28 
L2 5.40+1.18 5.76+1.35 0.14 
L3 5.1141.35 5.29 +1.53 0.31 
L4 7.53+1.76 7.7742.06 0.32 
Scaled by TBV 
Length 0.59+0.05 0.66 +£0.06 <10~4 
Width 0.21+0.02 0.22 +0.03 0.35 
Ll 0.053 +0.011 0.056+0.010 0.14 
L2 0.048+0.011 0.053+0.011 0.07 
L3 0.045 +0.012 0.048+0.013 0.20 
L4 0.067 +0.016 0.071+0.017 0.22 

Anterior 

Posterior 
(a) Bottom Top 

(b) 4 7 


Figure 7. (a) Overlaid average structures (blue: controls, red: patients) and (b) distance 
map between the two averages. 


3.3. Group mean difference map 


Figure 7 gives a descriptive visualization of the group difference. The overlaid average CC structures 
are shown in Figure 7(a). Figure 7(b) shows the signed distance map between the two average 
structures rendered on an overall average shape, where the negative distances indicate that the 
patients’ structure is outside the controls and the reverse is for the positive distances. The anterior 
is more projected and the posterior is more inward in the patients’ structure. The posterior body 
and isthmus are thinner in the patients’ structure. 


3.4. Statistical testing results 


The significance maps of the raw and corrected p-values are shown as color coded p-values in 
Figure 8. Smaller p-values indicate larger statistical significance. A two-tailed alpha level of 0.05 
is chosen as the significance threshold. Figure 8(a) is the significance map of raw p-values, which 
reveals significant difference in the anterior most, anterior body, isthmus and posterior bottom. 
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Top 
Bottom 
(a) Anterior Posterior (b) 
p=0 p>=0.05 


Figure 8. Significance map: (a) raw p-values and (b) corrected p-values. 


This result is an optimistic estimation. The permutation result is shown in Figure 8(b), which 
retains part of the significance in the raw map. The global shape analysis has a p-value of 0.008. 


4. CONCLUSION AND DISCUSSION 


Based on the new definition of autism subgroups, we investigate the CC abnormalities in essential 
autism, and exclude the complex autism group in order to find more homogeneous results. We 
use a semiautomatic method to segment 2D CC boundaries, and 3D surfaces are reconstructed 
by contour stitching. Regional volumes and additional phenotypic traits based on the oriented 
bounding rectangle of the CC are compared between patients and controls. A newly developed 
planar conformal mapping is used to parameterize the 3D surfaces. Shape comparison is conducted 
at each surface location using Hotelling T? test followed by a permutation test. 

In the volume comparison, no significance is found in TBV or ICV between patients and 
controls, but the non-significant trend for an increase in the brain volumes in patients is consistent 
with the most previous literatures [2, 3, 20]. The significant reduction in the total CC volume in 
the patients is consistent with [1]. We find significant reduction in each sub-region of the CC in 
the patients, but in [1] only a significant reduction in the anterior third was found. Besides the 
traditional volume test, we also conduct some tests based on the oriented bounding rectangle, 
which has not been done in the previous work. The length, width, and aspect ratio of the bounding 
rectangle are measured for comparison. We find significant reduction in the CC length and ratio 
in the patients, but the difference in the width is far from significant. This gives us some insight 
that the decrease in the CC volume is caused by the decrease in the anterior—posterior length more 
than the top—bottom length. 

In the shape analysis, the distance map reveals more projected anterior most and inward posterior 
bottom, as well as the thinner body of the CC in the patients. The raw significance map suggests 
significant differences in the anterior most, anterior body, isthmus and posterior bottom, while the 
corrected significance map retains part of the significance in each region shown in the raw map. 
Table II gives a summary of the findings in [1,6] and our study regarding the shape differences 
between autism and control groups. First, our results are similar to [1] in the posterior part, but 
the result in [6] is opposite. This may be accidental because there is no significance in their result. 
Second, our result on the anterior part is consistent with [6] but contrary with [1]. This may be 
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Table III. Comparison of the CC shape analysis in different studies. 


Body Anterior Posterior 
CC regions Shape Significant Shape Significant Shape Significant 
{1] Less arching No More inward Yes More inward Yes 
[6] Less arching Yes More projected No More projected No 
Our study Posterior thin Yes More projected Yes More inward Yes 


explained by the population because [6] and our study exclude both the complex autism while [1] 
does not. Both [1, 6] have found less arching in the body of the CC in autistic patients, but we 
have found the body thinner in the posterior half. He et al. [6] and our study both find significance 
in the shape changes of the body, although the changes are different. 

The results of this study need to be interpreted cautiously because of several limitations. First 
of all, the age range of our sample is very wide, which may have some effect on the CC measure- 
ments. Second, we exclude the complex autism group as it was done in [6]. However, with the 
improvements in the surface matching and statistical testing methods compared with [6], our results 
are more convincing than [6]. 

In conclusion, patients have a trend for a reduction in the volumes of the total CC and each sub- 
region. Further experiment with larger sample size needs to be done to confirm the significance of 
these size differences. New phenotypic traits based on the oriented bounding box are studied. The 
decrease in the anterior—posterior length is found to be the main reason for the CC size reduction. 
A new surface matching technique is used for shape comparison. The permutation results show 
significant difference in the anterior most, anterior body, isthmus and posterior bottom. Further 
studies on different subgroups within autistic patients need to be done, such as male vs female 
and complex vs essential. These findings will better explain the inconsistent results caused by the 
population. 
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